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Abstract 

An important aspect of real ferromagnetic particles is the demagnetizing field 
resulting from magnetostatic dipole-dipole interaction, which causes large par- 
ticles to break up into domains. Sufficiently small particles, however, remain 
single-domain in equilibrium. This makes such small particles of particular 
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interest as materials for high-density magnetic recording media. In this pa- 
per we use analytic arguments and Monte Carlo simulations to study the 
effect of the demagnetizing field on the dynamics of magnetization switch- 
ing in two-dimensional, single-domain, kinetic Ising systems. For systems 
in the "Stochastic Region," where magnetization switching is on average ef- 
fected by the nucleation and growth of fewer than two well-defined critical 
droplets, the simulation results can be explained by the dynamics of a sim- 
ple model in which the free energy is a function only of magnetization. In 
the "Multi-Droplet Region," a generalization of Avrami's Law involving a 
magnetization-dependent effective magnetic field gives good agreement with 
our simulations. 

FSU-SCRI-95-114 

PACS Number(s): 75.60.-d, 75.40.Mg, 05.50.+q, 75.10.Hk 
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I. INTRODUCTION 



The ability of single-domain ferromagnets to preserve an accurate record of past mag- 
netic fields has several important applications. Fine grains in lava flows preserve a record 
of the direction of the geomagnetic field at the time they cooled, giving valuable insight 
into continental drift and the dynamics of the earth's core.0 Of more direct technological 
importance is the potential application of single-domain ferromagnets to magnetic recording 
media, such as magnetic tapes and disks. 

During the magnetic recording process, different regions of the medium are briefly ex- 
posed to strong magnetic fields, so that each grain is magnetized in the desired directionJI 
Since each grain can in principle store one bit of data, a greater storage density could ideally 
be achieved by a medium containing many small grains than by one containing a few large 
grains. However, in order to serve as reliable storage devices, the grains must be capable of 
retaining their magnetizations for long periods of time in weaker, arbitrarily oriented ambi- 
ent magnetic fields — i.e., they must have a high coercivity and a large remanence. Since 
experiments show the existence of a particle size at which the coercivity is maximum (see, 
e.g., Ref. |3]), there is a tradeoff between high storage capacity and long-term data integrity 
which must give rise to an optimum choice of grain size for any given material. During both 
recording and storage, the relationships between the magnetic field, the size of the grain, 
and the lifetime of the magnetization opposed to the applied magnetic field are therefore of 
great technological interest. 

Fine ferromagnetic grains have been studied for many years, but until recently such 
particles could be studied experimentally only in powders (see, e.g., Ref. |]). This made 
it difficult to differentiate the statistical properties of single-grain switching from effects 
resulting from distributions in particle sizes, compositions, and local environments, or from 
interactions between grains. Techniques such as magnetic force microscopy (MFM) (see, 
e.g., Refs. and Lorentz microscopy (see, e.g., Ref. [D now provide means for overcoming 
the difficulties in resolving the magnetic properties of individual single-domain particles. 

The standard theory of magnetization reversal in single-domain ferromagnets is due to 
NeelS and Brown.EJ In order to avoid an energy barrier due to exchange interactions be- 
tween atomic moments with unlike orientations, Neel-Brown theory assumes uniform rota- 
tion of all the atomic moments in the system. The remaining barrier is caused by magnetic 
anisotropy,0 which may be either intrinsic or shape-induced. Anisotropy makes it ener- 
getically favorable for each atomic moment to be aligned along one or more "easy" axes. 
Buckling, fanning, and curling are, like uniform rotation, theoretical relaxation processes 
with few degrees of freedom and global dynamics.i'@ 

However, for highly anisotropic materials there exists an alternative mode of relaxation 
with a typically much shorter lifetime. Small regions of the phase in which the magnetization 
is parallel to the applied magnetic field (the "stable" phase) are continually created and de- 
stroyed by thermal fluctuations within the phase in which the magnetization is antiparallel 
to the field (the "metastable" phase). As long as such a region (henceforth referred to as a 
"droplet") is sufficiently small, the short-ranged exchange interaction with the surrounding 
metastable phase imposes a net free-energy penalty, and the droplet will, with high proba- 
bility, shrink and vanish. Should the droplet become larger than a critical size, however, this 
penalty will be less than the benefit obtained from orienting parallel to the magnetic field, 
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and the droplet will with a high probability grow further, eventually consuming the grain. 
The nature of the metastable decay thus depends on the relative sizes of the grain, the crit- 
ical droplet, the average distance between droplets, and the lattice constant, as discussed in 
detail, e.g., in Refs. |14|-[T7|. For systems in which short-range interactions dominate, Fig. [I] 
sketches the four regions of the space of magnetic fields and particle sizes distinguished by 
different behaviors during metastable decay. For a more complete, recent review of droplet 
theory, see Ref. [T7|. 

It is important to understand the difference between a droplet and a domain.ll! Although 
they are both spatially contiguous regions of uniform magnetization, a domain is an equilib- 
rium feature whereas a droplet is a strictly non- equilibrium entity. The domain structure of 
two-dimensional dipole systems has been extensively investigated-EHnl The magnetostatic 
dipole-dipole interaction produces a demagnetizing field, which results in the stabilization 
of a domain structure in large ferromagnetic particles at equilibrium. 

The purpose of the present paper is to study the effects of long-range dipole-dipole 
interactions on the nonequilibrium phenomenon of magnetization switching in single-domain 
ferromagnetic particles. Towards this end we employ a simplified model in which particles 
in equilibrium can have only one or two domains, and we emphasize the single-domain case. 

Detailed descriptions of both the static and dynamic properties of fine ferromagnetic 
grains have typically been formulated from micromagnetic studies.il This method involves 
coarse-graining the physical lattice onto a computational lattice and then solving the partial 
differential equations for the evolution of magnetic structures on the computational lattice. 
Although micromagnetics provides a good treatment for the anisotropy and demagnetizing 
fields, it treats thermal effects rather crudely, usually just by making the domain- wall en- 
ergy temperature-dependent. A somewhat better approximation for thermal fluctuations 
within the underlying differential equations is to include small fluctuations using a Langevin 
noise term.!§ A better treatment for thermal and time-dependent effects is therefore Monte 
Carlo simulation (see, e.g., Refs. P?|-[ZI5|). Even when the physical phenomena can be accu- 
rately simulated, however, it will be difficult to understand the results without an adequate 
theoretical basis. 

Because of its simplicity, the kinetic nearest-neighbor Ising model has been extensively 
studied as a prototype for metastable dynamics (see Ref. [17| and references cited therein). In 
particular, square- and cubic-lattice Ising systems with periodic boundary conditions have 
been used to study grain-size effects in ferroelectric switching.0'E3 A related one-dimensional 
model has been used to study magnetization reversal in elongated ferromagnetic particles.il 
In a recent article,^ we applied statistical-mechanical droplet theory and Monte Carlo sim- 
ulations of two-dimensional Ising systems to obtain a qualitative approximation for the 
dynamical behavior of real single-domain particles magnetized opposite to an applied field. 
In so doing, we made several simplifying approximations, one of which was the absence of a 
demagnetizing field. In this article we consider the effect of a small demagnetizing field on 
Ising systems, and compare our analytic calculations with simulations of two-dimensional 
Ising systems. Specifically, for systems in the Stochastic region (discussed in Sec. [TIJ), the 
demagnetizing fields we consider must be sufficiently small that the system is consists of a 
single domain in equilibrium, whereas in the Multi-Droplet region (discussed in Sec. fV\ ) it 
is sufficient to have for demagnetizing field to be much smaller than the applied field. Some 
preliminary results were presented in Ref. 
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The organization of this paper is as follows. In Sec. |H| we define the model and numerical 
methods employed in this paper. In Sec. Ill we discuss the Stochastic region in terms of an 
approximate free-energy functional and give some numerical results. In Sec. [IV] we generalize 
Avrami's Law ,HHEI which describes magnetization switching in the Multi-Droplet region, to 
include the effects of the demagnetizing field, and we compare the analytical results to 
numerical simulations. Section [Vj contains conclusions and discussions. 



II. MODEL AND NUMERICAL METHODS 



The standard Ising model is defined by the Hamiltonian 

Ho = — J E s i s j ~ HL d m , 



(1) 



where Sj = ±1 is the z-component of the magnetization of the atom (spin) at site i, J > is the 
ferromagnetic exchange interaction, and H is the applied magnetic field times the single- 
spin magnetic moment. The sum Yluj) runs over all nearest- neighbor pairs on a square 
(generally d- dimensional hypercubic) lattice of side L. In this work we do not study the 
effects of grain boundaries, so periodic boundary conditions are imposed. The dimensionless 
system magnetization is given by 



m 



(2) 



where the sum is over all L d sites. The lattice constant is set to unity. 

Addition of dipole-dipole interactions gives a total Hamiltonian (SI units) 



n 



V 



H + 



47T 



E 



SiSj 



(3) 



where M is the saturation magnetic dipole moment density and is the vector from 
site i to site j. Unfortunately, however, the last sum in Eq. @ slows down Monte Carlo 
simulations significantly, which is problematic if a large number of realizations are desired 
for good statistics, as is the case in nonequilibrium studies. The last sum also would make a 
perturbative expansion in the demagnetizing field (adjustable by changing M or the sample 
shape) difficult. We therefore instead use the simpler Hamiltonian 



n D = H + L d Dm 2 



(4) 



where D is a function of the crystal symmetry, the shape of the system, and M. Equations ([3]) 
and (|4]) are equivalent for general ellipsoids uniformly magnetized along a principal axis. 
For the special case of a perpendicularly magnetized plane with square-lattice symmetry, 
D = ^fi M 2 . For non- uniformly magnetized systems, Eq. (Q) amounts to a mean-field 
treatment of the effects of the dipole-dipole interactions. 

For systems with periodic boundary conditions, the exchange and dipole terms of Eq. (|J) 
are equal when the system size is given byEl 
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where cx^T) is the surface tension along a primitive lattice vector in the limit L — > oo and 
m sp (T) is the spontaneous magnetization. For the two-dimensional Ising model, cr^T)!! 
and m sn are known exactly. The length scale on which we would expect a transition 
from a single-domain to a multi-domain equilibrium structure is approximately Lp. 

The selection of the Ising model is equivalent to requiring a very large (infinite, in 
fact) anisotropy constant. Although magnetic materials used in magnetic recording media 
require comparatively large anisotropy constants,! the microscopic anisotropy tends to be 
much smaller than the exchange energy. However, the role of the anisotropy is enhanced 
by coarse-graining. Simplicity is our main reason for choosing the two-dimensional Ising 
model with periodic boundary conditions, particularly since manv equilibrium properties of 
the two-dimensional Ising model in zero field are known exactlyBa and since the kinetics 
of metastable decay has been extensively studied for this model.0 As a result, our model 
systems more closely resemble ultrathin magnetic films with perpendicular magnetization 
than magnetic grains. A more realistic simulation of three-dimensional grains is planned 
for later study, but we emphasize that we expect that droplet theory applies to almost any 
spin model with high anisotropy. Accordingly, equations are written in forms appropriate 
for arbitrary dimensionality d, even though simulations are only carried out for d = 2. 

The relaxation kinetics is simulated by the single-spin-flip Metropolis dynamic with up- 
dates at randomly chosen sites. A rigorous derivation from microscopic quantum Hamiltoni- 
ans of the stochastic Glauber dynamic used in Monte Carlo simulations of Ising models has 
been established in the thermodynamic limit ;0 both the Glauber and Metropolis algorithms 
are spatially local dynamics with non-conserved order parameter (the dynamic universality 
class of Model A in the classification scheme of Hohenberg and Halperincll) and are therefore 
expected to differ only in non-universal features. The Metropolis dynamic is realized by the 
original Metropolis algorithm^ and the n-fold way algorithm.^ (For a discussion on the 
equivalence of the dynamics of these algorithms, see Ref. The acceptance probability in 
the Metropolis algorithm for a proposed flip of the spin at site a from s a to — s a is defined 
as W(s Q — > —s a ) =min[l, exp(— 0AE a )], where AE a is the energy change due to the flip and 
= ksT is the temperature in units of energy. The n-fold way algorithm is similar, but 
involves the tabulation of energy classes. First an energy class is chosen randomly with the 
appropriately weighted probability. A single site is then chosen from within that class with 
uniform probability and flipped with probability one. The number of Metropolis algorithm 
steps which would be required to achieve this change is chosen from a geometric probability 
distribution,^ and the time, measured in Monte Carlo steps per spin (MCSS), is incremented 
accordingly. The n-fold way algorithm is more efficient than the Metropolis algorithm at 
low temperatures, where the Metropolis algorithm requires many attempts before a change 
is made. 

Because we are using single-spin-flip dynamics, the magnetization can only change by a 
small amount from one time step to the next. The dynamical effects of the demagnetizing 
field thus depend only on the change in the magnetic part of Hp between adjacent values 
of the magnetization. In this way it is possible to define an effective magnetic field 

d 

H cS (H, D,m) = — [Hm - Dm 2 
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= H - 2Dm 



(6) 



The effective magnetic field is thus site-independent. This fact makes analytic considerations 
significantly easier and is our principal reason for using Eq. (Q) rather than Eq. (^) as the 
Hamiltonian. 

We study the relaxation of the dimensionless system magnetization starting from an 
initial state magnetized opposite to the applied field [rn(t — 0) — +l,H<0]. This approach 
has often been used in previous studies, e.g. in Refs. |T6| , f45| . For the temperatures employed 
in this study, the equilibrium spontaneous magnetizations in zero field are close to unity, 
with 0.95 < m sp < 1. Since the applied field is negative (and generally small), the stable 
magnetization is approximately m st « — m sp and the metastable magnetization is m ms m 
+m sp . We use as an operational definition of the lifetime r of the metastable phase the 
mean first-passage time to a cutoff magnetization m = 0: 

r=(t(m=0)>. (7) 

It has been observed^ that the qualitative results discussed below are not sensitive to the 
cutoff magnetization as long as it is sufficiently less than m sp . Our choice of m = as the 
cutoff facilitates comparison with MFM experiments, which are only capable of measuring 
the sign of the particle magnetization. 

In this paper a numerical subscript indicates the coefficient in a Taylor expansion in 
D. For example, a quantity X may be expanded X = X + X\D + X 2 D 2 + . . .. There are 
three exceptions to this rule: (1) The subscripts in Eq. (|24f) refer to an iterative process for 
evaluating a continued fraction. (2) The subscripts on S (T) and Si(T) [Eq. fl2"E|)] indicate 
an expansion in H 2 and are kept for continuity with Ref. |33]. (3) Dummy variables in the 
Appendix [e.g., Eq. ( |A2|) 1 may have numerical indices as a matter of convenience. 



III. THE STOCHASTIC REGION 

It has been showniHil that the dynamics of metastable decay in the standard two- 
dimensional Ising model for sufficiently weak applied field can be semiquantitatively de- 
scribed by a mean-field-like dynamic in which the free energy is a function only of the system 
magnetization. Under these circumstances switching is abrupt, with a negligible amount of 
time being spent in configurations with magnetizations significantly different from m ms or 
m st . Switching is also a Poisson process, with the lifetime of a metastable phase given by 
the typical Arrhenius form 

rocexp(/3AF) , (8) 

where AF is the free-energy barrier that must be crossed in the decay process, or by a simple 
generalization of Eq. (||) if more than one equivalent decay path is present [see Eq. (^5|)]. 
This phenomenon, in which the entire system behaves as though it were a single magnetic 
moment, is known as superparamagnetism.00 The standard deviation of the switching time 
for an individual grain is approximately equal to the mean switching time, r. Because of the 
random nature of switching in this region, it has been calledSEl the "Stochastic" region. 
The Stochastic region is the union of the "Coexistence" region and the "Single-Droplet" 
region (discussed below). 
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In the spirit of Refs. |46|-|48| we construct an approximate restricted free-energy function 
F(m) for the entire system and use Eq. (||) to illustrate the H- and D-dependence of the 
lifetime: 



F(m) = L d Dm 2 + min{F u , + (m), F u ,_(m), F d ,+(m), F d ,_(m), F sl (m)} - F sl (m = 0) , 



(9) 



where F s i(m) is the free energy of a system composed of two "slabs" with magnetizations 
near ±m sp [Fig. |](a) illustrates a "slab" configuration], F di ±(m) is the free energy of a system 
with a single droplet with magnetization near =Fm sp in a background with magnetization 
near ±m sp [Fig. |](b) illustrates a "single droplet" configuration], and F Uj ±(m) is the free 
energy of a system in a "uniform" phase near m = ±m sp . Figure |3] illustrates f3F(m) for 
L<^:L D , L = L D , and L^>L D . 

We approximate the free energy of a system in a "uniform" phase by 

F Ui± (m) = -L d iYm + L^x^m T ™s P ) 2 , (10) 

where \ is the equilibrium susceptibility per spin. Since an exact solution for the two- 
dimensional Ising model in a magnetic field has not yet been found, we use instead an 
estimate from a series expansion,B2l so that for T = 0.8T C , x^0.05J _1 . 

Minimizing F u _{m) + L d Dm 2 yields the stable magnetization m st , which is the location 
of the global minimum of F(m) for LkL^: 

-m sp + Hx M 1 v 

m st w 11a 

1 + 2D X 

(remember, H <0). Likewise, for L<^Ld the next-lowest minimum of F(m) is obtained my 
minimizing F Uj+ (m) + L d /}m 2 : 

m sp + Hx /-, -, i \ 

771ms ~ ~ • Hb 

1 + 2D X 



Equation ( 11a ) is valid for a wider range of H than is Eq. ([Llb|) . We shall refer to m ms as 



the "metastable magnetization" and its basin of attraction as the "metastable phase" , since 
for systems of interest (L<Ld) the length of time required for a system initially prepared in 
the metastable phase to escape to the stable phase is much longer than any other timescale. 
Note, however, that other, shorter-lived metastable phases may exist (discussed below). 

In cases where the magnetization differs significantly from m sp (— m sp ), a lower free 
energy can often be obtained by segregating the system into a single localized "droplet" 
with magnetization near m s t (m ms ) in a background with magnetization near m ms (m s t)@ 
Specifically, the droplet free energy is approximated by 



F dt+ {m) = Q 

and 

Fa Jm) = Q 



da^Rl' 1 + (m ms -m st )HR d + - L d Hm ms (12a) 



da^Rt 1 - (m ms - m st )HR d _] - L d Hm st (12b) 



subject to m ras >m>m st . Here 
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R, = Q-V d L 



m m « - m 



l/d 



•■ ^ms ^st 

is the radius of a droplet of "down" (stable) spins in an "up" (metastable) background 

l/d 



(13a) 



r_ = fr 1/d L 



m — m st 



(13b) 



is the radius of a droplet of "up" (metastable) spins in a "down" (stable) background. 
For T 96 0, the droplet shape can be found from a Wulff construction. The quantity Q, 
which gives the volume of the droplet via V = QR d , can be found to arbitrary precision for 
the two-dimensional Ising model by numerically integrating over the exactly known surface 
tension. HH 

Lastly, near m = the circumference of the droplet becomes larger than twice the cross- 
section of the system, and the lowest free energy is obtained by segregating the system into 
two slab-like configurations." The corresponding slab free energy is approximated by 



Comparison with Eq. 



F A (m) = 2L d - 1 (T 00 - L d Hm. 
]) shows that F(m)=F s \(m) for md s ,+ > m > m^s 



(14) 



whercE10 



™ds 



m st )Q-^ (!) > 



and 



//Ms._ = met + (mms - m^n- 1 ^^ (-) < 



(15a) 



(15b) 



Note that for D > H / ^2m<is,-) , a local minimum of F(m) occurs for a slab configuration 
at m = H/ (2D). For L < Ld it is a metastable phase, but for L>Ld and D > H/m st it is 
the global minimum of F(m) and hence the true stable phase (see Fig. |^). Other interesting 
features can be obtained by solving (d/dm) L d Dm 2 + F^±(m) = 0, which can in general 
be done only numerically. This reveals that near L = Ld short-lived metastable phases can 
exist for m>md s ,+ or m<md s ,-- 

For L < Ld, the system enjoys true coexistence at zero applied field between two 
degenerate equilibrium phases with magnetizations m ms and m st . This leads to the 
identification^!'!!! of a "Coexistence" (CE) region within which F(m ms ) mF(m st ). Within 
the CE region, the free energy barrier for tunnelling from the metastable phase to the stable 
phase is approximately the same as the free energy barrier for tunnelling from the sta- 
ble phase to the metastable phase, so the decay process is both stochastic and reversible. 
Specifically, for L <^ L D , the life^me of the metastable phase is given by Eq. (0) with 
AF = F(m ds , + 



F(m ms ), so thatl 



t(L,H,T) w A(T) exp{(3 



2a 00 (T)L d - 1 - L d \H\(m ms - m ds , + ; 



L d D(ml, 



m 



ds,H 



(16) 
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where A(T) is a non-universal prefactor. 

For Lz^Ld, the maximum of F(m) occurs not at mj s>+ but at a larger magnetization 
corresponding to a single critical droplet. The size of this droplet, however, is strongly 
dependent on the system size L. This part of the Coexistence region is further complicated 
by the increasing importance of the metastable phase at m = H/ (2D) and the aforementioned 
possibility of metastable phases near m = md s ,±- 

For H^>D the maximum of F(m) may correspond to a critical droplet the size of which 
is nearly independent of system size since it is determined by the applied field rather than 
the demagnetizing field. The applied field at which the CE region crosses over into this 
"Single Droplet" (SD) regionElS is called0'@ the "Thermodynamic Spinodal" (-ffrhSp)- A 
useful estimate for this crossover is (d/dm) L d Dm 2 + F dj+ (m) =0, yielding 



Ms, 



\H ThSp \ « L-^W-D ( f]^ _ 2Dm dSt+ . (17) 



m ms ~ m st 



A slightly different estimate was used in Ref. |Tj and Ref. |33], but the two estimates are 
approximately equal. 

In the SD region the first critical droplet to nucleate almost always grows to fill the 
system before any other droplet nucleates. The average time required to nucleate the first 
droplet can be estimated from Eq. where the free-energy barrier is determined from 
Eq. ( [T0|) and Eq. ( 12a ). Since the SD region is also the region of weak H and D, we can 



obtain a good approximation by neglecting terms of 0(H e gx)- Then the magnetization in 
the metastable background is m ms = m sp , and inside the droplet it is m st = —m sp . In terms 
of the droplet radius R, the difference between the free energy of a system containing one 
droplet and that of a uniform metastable system can then be written as 

AF(R) = dtta 00 (T)R d - 1 - 2m sp (\H\ + 2Dm sp ) ttR d + ADm 2 sp L- d n 2 R 2d . (18) 

Differentiating with respect to R, we find the implicit equation satisfied by the critical 
droplet radius: 

R C (T, H, D) = (rf ~ lKo(T) N| , (19a) 



where 



\H eStC (H,D)\ = \H\ + 2Dm sp \l - 2VL (R c /L) d ] (19b) 



is the effective field evaluated at the magnetization of a system containing a single, critical 
droplet. Note that H e s c depends on \H\ and D explicitly, as well as implicitly through 
R c . For D = 0, Eq. flTjp reduces to the standard expression for RM By inserting R c into 
Eq. flTSP one finds that the free-energy barrier corresponding to the critical droplet is also 
simply given by a standard expression Jll in which \H\ has been replaced by |if e ff,c| : 

^ - HOT), ■ < 2 °> 
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where0 



S (T)=/?fi [^(T)]' 



'd 



d-l 



2m. 



(21) 



sp 



Note that S (T) is completely defined byauantities that for the two-dimensional Ising model 
are either known exactly (a^ and m sp )H'ii or can be obtained by numerical integration of 
exactly known quantities (fi).HI! 

The above results indicate that in order to obtain the nucleation rate for nonzero D, one 
only needs to determine |-ffeff,ci- This can easily be done to arbitrary numerical precision via 
a rapidly convergent, generalized continued- fraction expansion as follows. 

Let x = 2Dm sp /\H\ be the reduced demagnetization field and V(x) = 2Q(R C / L) d be the 
volume fraction occupied by the critical droplet. Then 

|# effjC | = \H\ {1 + x [1 - V(x)}} = \H\ y(x) , 

and V(x) is given by the generalized continued-fraction expansion, 

Vo 



(22) 



V{x) 



( 



1 + x 



1 - T 



V 



\ 



1+cc 1- 



(23) 



where Vo = V(0). This expansion can be evaluated to desired precision by the recursion 
relation, 



Vn 



[l + x(l- K-i)] ' 



(24) 



(Here the subscript is proportional to the order of a rational-function approximation to the 
generalized continued fraction, rather than denoting the order of a term in a power series in 
D, as elsewhere in this paper.) 

The lifetime in the SD region can be given in terms of the nucleation rate per unit volume 
Jby0 



L d I 



-i 



For D = and x~0, J has been shown to be given by0 

I(T,H)^B(T)\H\ K exp{-\H\ 1 - d [E (T) + E 1 (T)H 2 ]} , 



(25) 



(26) 



where B(T) is a non- universal prefactor, and K is believed to be 3 for the two-dimensional 
Ising model and —1/3 for the three-dimensional Ising model.B (The subscripts on S (T) and 
Si(T) indicate an expansion in H 2 rather than in D, and are kept for consistency with the 
notation in Ref. |33|.) The quantity So(T) is given by Eq. (21), and we determine Si(T) from 
a numerical fit to the if-dependence of the lifetime. Note that the lifetime obtained from 
Eqs. (|25| ) and (|26|) is quite similar to the Arrhenius form with the free-energy barrier given 
by Eq. (|20|) . The only differences are the prefactor \H\~ K and the term Ei(T)H 2 in the 
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exponential, which are due to surface fluctuations on the droplet and to higher-order terms 
in a field-theoretical calculation of the free-energy barrier, respectively.0 We generalize to 
the D^O case by assuming that the nucleation rate in the SD region is given by Eq. ( |26|) 
with | if | replaced by \H c{ic (H, D)\, as we have already shown in Eq. (|20|) for the dominant 
term in the free-energy barrier, AFsd- The resulting expression for the relative lifetime for 
nonzero D is then given by 



T ( x ) 
r(0) 



y(x) 



-K 



exp {-S |if | 1_ 



l-y(x) 



l-d 



-EAH 



3-d 



i-y(x) 3 - d }} 



(27) 



where y(x) is defined in Eq. 



This result is shown in Fig. [|, together with MC data, 



for d=2, T=0.8T C , H=0.2J, and L=10. Except for Si, the parameters needed to evaluate 
t(x)/t(0) are known exactly or numerically exactly for the two-dimensional Ising model: 
^o=0.2399(l) and S =0.5062(l). The value of Hi used in Fig. |, Si = 9.1(3) J' 1 , was 
obtained in Ref. |33] from the \H\ dependence of r for D = in the SD region. Thus the 
good agreement seen in Fig. |4] between the simulation data and the theoretical prediction is 
not the result of a fit to the data, but is determined entirely from quantities measured with 
D = 0. 



IV. THE MULTI-DROPLET REGION 

For sufficiently strong fields or large systems, decay occurs through many weakly interact- 
ing droplets in the manner described by Kolmogorov,H Johnson and Mehl,@ and Avramilll 
Such decay is "deterministic" in the sense that the standard deviation of the switching time 
is much less than its mean (see Ref. [!3| for details). The crossover between the SD region 
and the Mult i- Droplet (MD) region has been calledffl the "Dynamic Spinodal" (DSp). 
Since the standard deviation of the lifetime is equal to its mean in the Stochastic region, we 
estimate this crossover by the field ifi/2 at which 

y/(t*(m = 0))-T* = ^. (28) 

For asymptotically large L, if DSp ~ (1/ lnL) 1 ''^ -USE however, prohibitively large system 
sizes may be required before this scaling form is observed. At a sufficiently high field, 
nucleation becomes much faster than growth and the droplet picture breaks down. The 
crossover to this "Strong-Field" region (SF) has been the "Mean-Field Spinodal" 

(MFSp). A conservative estimate for this crossover field is obtained by setting 2R C =1: 

| ffMFSp , „ ('-D-MD , (29) 

Little is understood quantitatively about the SF region. 

In this section we generalize Avrami's LawETLlI to systems with nonzero D. Avrami's 
Law gives the volume fraction of the metastable phase (or equivalently, the magnetization) 
for systems in which droplets nucleate with a constant rate (per unit volume) Iq and grow 
at constant velocity Vq without interacting except for overlaps. The generalization we make 
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below is for a nucleation rate and velocity that depend on the magnetization, and through 
it on time. 

The time-dependent mean system magnetization m(t) is given by 



m(t) 



m st )e * w + m st 



Here 



$(t) = /' ' I(t')V(t',t)dt' 
Jo 



is the mean volume fraction of droplets (uncorrected for overlap) and 

i d 



v(h,t 2 ) = n 



>2 



v{t)dt 



(30) 



(31) 



(32) 



is the volume occupied by a droplet which nucleates at time t\ and grows with a 
time-dependent radial velocity v(t) until time t 2 - Here v(t) is the (nonuniversal) 
temperature-dependent radial growth velocity of a droplet, which under an Allen-Calm 



approximatio 



is proportional to the effective field in the limit of large droplets: 



v{t) H eS (H,D,m{t)) 



(33) 



The time-dependent nucleation rate is given by I(t) = I[T,H eS (H,D,m(t))] from Eq. 
Note how this differs from the D-dependent nucleation rate in the SD region. In the SD 
region, the .D-dependence of the nucleation rate comes from the change in system magne- 
tization from the nucleation of a single critical droplet. In the MD region, by contrast, we 
ignore the change in system magnetization due to the nucleation of a single droplet (since 
L 3> R c ), and the -D-dependence of the nucleation rate comes from the change in system 
magnetization due to an ensemble of droplets. 
For D = 0, Eq. (pCf) becomes 



m (t) 



ra s t,o)e * oW + m sti0 , 



where 



$o(t) = In 



m ms,0 — m st,0 

|m st ,o| 



TO J 



(34) 



(35) 



so that To is the first-passage time to m = [Eq. (0)]. Specifically, r is given 



TO 



d+1 



where 



(d + l)\nzQ 



\m s t,o\ 



2 . 



(36) 



(37) 



In Eq. ( p^) and elsewhere in this section, the estimate for the metastable magnetization 
given by Eq. (|llb| ) does not suffice. Ramos et al. have estimated m ms (H) by extrapolating 
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m{t) back to t = assuming that m{t) is correctly described by Avrami's Law.EI (Since the 
initial condition is m = 1 rather than m = m ms , the earliest times must be discarded.) 
These estimates are shown in Fig. [5]. We fit a smooth curve through the data, insisting that 



sp 



and (d/ l dH)m ms \ H=0 = x- The smooth curve allows us to estimate the change 



in m ms due to D, but the correct form of m ms (H) is not known. 

For D = 0, the standard deviation of the time-dependent magnetization has been shown to 
vanish with increasing system size as L _d//2 ,ilS a feature which is shared by the more general 
case of D > 0. In fact, realizations in which the magnetization chances to decay more rapidly 
than average will experience a weaker effective field [Eq. (H)], and realizations in which the 
magnetization decays more slowly than average will experience a stronger effective field. 
These effects combine to cause systems with D > to have even smaller standard deviations 
in their time-dependent magnetizations than corresponding systems with D = 0. 

To first order in D, the effective magnetic field [from Eq. @] is given by 



H eS (H,D,m{t))^H-2Dm {t) 



(38) 



since any D-dependent terms in m(t) will lead to only higher-order corrections. We will 
expand $(t) to first order in D, so that we can use the known value of mo(t) instead of the 
unknown value m(t) on the right-hand side of Eq. (|3TD . In order to perform the expansion 
correctly, we must expand I(t) and V(ti,t 2 ) to first order in D. Specifically, the total volume 
fraction (uncorrected for overlap) of droplets after a time t is given by 



where 



$(t)«$ (*) + $i(t)£> 



(39) 



(40) 



and 



Jo 



/ h{t')V {t' ^t' . 



(41) 



(42) 



Straightforward but cumbersome mathematics shows (see the Appendix) 



Mt) « 2(d)|F|- 1 Kt ) o|<| -%(t) + Z j [$ (t)] d/(d+1) 7 



1 



+ 



^eC; 1 )!-!)^^)]^ 



A 



d + r 
~-(k + i) l 



d+l d+1 



and 



«A(iy)|m st , |^oi:|(f) (-l) d - fc [$o(t)f /(d+1) 7 
- $o(t)l , 



d + l 



,$o(t) 



$o(*) 



(43) 



(44) 
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where 



X(H) 




(d-l)So + (d-3) S x # 5 



(45) 



comes from differentiating Eq. ( [ZED with respect to H, 7 denotes the incomplete gamma 
function 



and A is given by 



/■X 

7(0,2;) = / y^e^dy 

JO 



*4(a,6,a;) = / y a 1 {b,y)dy 
Jo 



(46) 



(47) 



Both 7 and A are easily evaluated by Taylor series. 

We can insert $y(t) and $j(i) into Eq. ((40]) to find $i(t) and then use Eq. (|30|) to 
evaluate 



mi(t) 



[171 



ms,l 



m st A e 



-*o(t) 



ms,0 



- mst ,„)e- $ »W$ 1 (t). 



(48) 



Finding the second-order terms in D proceeds along parallel lines. Although it is possible 
to find an analytic expression for rri2{t), this expression is tedious to derive and unenlight- 
ening. Furthermore, enough approximations have already been introduced to make the 
significance of an analytic expression for rri2(t) suspect. Consequently, we estimate m(t) 
by integrating Eq. fl3"T|) numerically with the effective field H — 2D[mo(t) + Drrii(t)]. The 
resulting estimate we denote mj(t), and it should be approximately correct to 0(D 2 ). If 

necessary, this process can be iterated to find successively better approximations for m(t). 

Figure ^ shows the time-dependence of the the magnetization, both as approximated 
above and as simulated by Monte Carlo. For Hi we have used Hi = 3.0(3), as determined in 
Ref. ^ from the \H\ dependence of r for D = in the MD region. Note that there is good 
agreement between the simulation results and the approximation after an initial relaxation 
into the metastable phase, and that the modification to m(t) resulting from higher-order 
terms is minor. 

In order to calculate the effect of D on the lifetime, we start with m(r) =0 and expand 
both m and t in D. Collecting terms and discarding all terms of higher order than D 2 , we 
find 



= {m (r )} + |mi(r ) + n 
dm 1 



dm, 



dt 



D 



+ { m 2 (T ) + n 



dt 



+ r 2 



TO 

dm. 



TO 



dt 



1 2 d 2 m 



TO 



V 1 dt 2 



D 



(49) 



TO - 



Since the quantities in the braces are independent of D and Eq. ( J4"9"D is true for all small D, 
by necessity mo(T ) = 0, 



Tl = -mi(To) 



dm 




-1 









(50) 
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and 



r 2 = ~ 



dm! 



dt 



TO 



+ 2« 



d 2 mo 






dm 




dt 2 


TO- 




dt 


no- 



(51) 



Equation ([50]) is readily evaluated because <&o(t) is a simple function. It is likewise simple 
to calculate (d 2 m /dt 2 )| T0 for use in Eq. fl5"T|). We have not actually solved for m 2 (t), but 
for small D we can use 



7712 (t) ~ D- 2 {mr(t) - [m (t) +Dmi(t)]} 



(52) 



Finally, (dmi/dt)| ro can be evaluated by differentiating Eq. (f|3|) and Eq. ([PJ) with respect 
to t: 



2(rf)|^|- 1 |m st 



(d+ l)lnz 



<=T0 



To 



X • 



-1 + - — , (mzo) /v ; 



1 



7 



ln^o 



2(1 



d-1 



D 2 S 



A 



f-(k + l) 



d- 1 



d + r 

d - k k(lnz f- d - 1)/id+1) 



d+l 'd+1 



lnz 



(53) 



and 



df 



t=T0 



A(if)|m st , 
x (lnz 



(d + 



o 

(fc-d-l)/(d+l) 



l)lnz / 



i + 5Tl£{(J)(-i)"* 



7 1 



d + l 



lnz 



(54) 



This yields (d$i/dt)| T0 from Eq. (|40|) . Once this is known, differentiating Eq. fl48f ) is trivial, 
and r 2 can easily be evaluated. 

Figure [7] shows r vs. D for two different values of H . The agreement between the 
theoretical curves and the Monte Carlo data is again excellent. Once again, the theoretical 
curves are not fits to the simulation data, but use only parameters determined for D = 0, 
namely, m ms (if)S and Hi(T)Ji 



V. DISCUSSION 

Due to the importance of magnetic recording technologies in modern society, magnetic 
relaxation has been a subject of study for many years. However, even the equilibrium 
thermodynamics of magnetic materials is very difficult to predict from first principles and 
generally has to be approximated from simpler models (see, e.g., Ref. |>|). As a result, 
the most popular method for theoretical investigation of magnetization reversal involves 
setting up and solving differential equations on a lattice obtained by course-graining over 
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the microscopic crystal lattice. This method, known as micromagneticsj! often gives very 
good results, particularly for equilibrium studies or for multi-domain particles. However, 
micromagnetic calculations take thermal effects into account only crudely. 

An alternative method is to treat the statistical mechanics carefully, making simplifi- 
cations to the model until it can be well understood. This is the approach we take. In 
Ref. |33| we showed that both the switching field and the probability that the magnetization 
is greater than zero, calculated from Monte Carlo simulations of the two-dimensional Ising 
model, are qualitatively similar to the same quantities measured in isolated, well character- 
ized single-domain ferromagnets by techniques such as MFM. Since statistical-mechanical 
droplet theory successfully explains the Ising model simulations, it is plausible that droplet 
theory could also be applied to the experimental particles. 

In this article we consider the effect of the magnetic dipole-dipole interaction, which was 
neglected in Ref. [33|. By treating the dipole-dipole interaction in a mean-field approximation, 
we are able to calculate droplet-theory predictions for the lifetime for systems in which 
magnetic decay occurs by means of a single droplet (Fig. f|). We also obtain both the time- 
dependent magnetization (Fig. |6]) and the lifetime (Fig. |7|) for systems in which magnetic 
decay occurs through the action of many droplets. In all of these figures, all parameters 
were determined by measurements at D = 0, so the excellent agreements are not the result 
of curve fitting. 

It should be pointed out that the droplet-theory predictions made in both the Single- 
Droplet and Multi-Droplet regions are large-droplet approximations. Since a critical droplet 
in the two-dimensional Ising model at T = 0.8T C and \H\ = 0.3 J consists of approximately six 
overturned spins, it is quite remarkable that these expressions give the good approximations 
they do. 

In Ref. [27], Kirby et al. use the two-dimensional Ising model with mean-field dipole-dipole 
interactions to simulate Dy/Fe ultrathin films which they have observed experimentally, 
obtaining good agreement. The analytic results from Sec. [TV] are directly applicable to such 
films, although different values of T, H, and D must be chosen to make the comparison. 
However, care should be used in applying the results. Specifically, if \Dm sp /H\ is not 
small, more numerical iterations of the type described in Sec. [TV] will be necessary, and if 
\H e g\ > |i?MFSp| [Eq. (|29"D 1, droplet theory may not be applicable. 

It is interesting to note that the mean-field demagnetizing field we have used does not 
change the qualitative behavior of the switching field as a function of system size from that 
which was studied in Ref. [J3] and sketched in Fig. [I|. Although the values of the switching 
field are reduced, a peak in the switching field still occurs near the thermodynamic spinodal, 
as can be seen by comparing Eq. (|H]) and Eq. (P5|). Furthermore, the switching field remains 
roughly independent of L in the MD region. What is noticeably absent is any feature at 
L = Lo'i the switching field shows features due to transitions in the dynamics, but not 
transitions in the statics. 

Even with the addition of the demagnetizing field, the Ising model remains too crude a 
model for magnetism to describe quantitatively real magnetic materials, except perhaps for 
some ultrathin films as noted above. Heterogenous nucleation at boundaries, quenched dis- 
order, and more realistic anisotropies are therefore planned as the subject of future studies. 
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APPENDIX: AVRAMI'S LAW FOR D>0 

In this Appendix we give some of the steps that have been omitted for clarity in Sec. [TV]. 
Beginning with the Allen-Cahn approximationEaEl for the radial growth velocity [Eq. ([33])] 
and with the effective magnetic field to 0(D) [Eq. (§)], we find 



v(t) « v 



2D 

W\ 



m (t) 



(Al) 



Substituting Eq. (|A1|) into Eq. (0), we find 




2D 

l + j^m (t) dt 
d\ 2D 



(h - h) d +^^R(t 2 - Uf- 1 f* mo(t)dt 



(A2) 



This enables us to make the identifications 



V (t u t 2 ) = Qv${t 2 - hY 



and 



Vifa, h) = 2nv d ^) li^- 1 (t 2 - hf- 1 J** m (t)dt . 
A Taylor expansion of the nucleation rate [Eq. (|26f) 1 

d\H cS \ 



(A3) 



(A4) 



I(t) ^I + 



d/ 
dD 



D = I + 



Io + 



dl 

d\H eS \ 

dl 
d\H> 



eff 



dD 



2m Q {t) 



D 



D 



(A5) 



yields 
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h(t) = X(H)I m (t) 



(A6) 



where A (if) is given by Eq. (p)|). 

The evaluation of and $y(t) can be followed more easily if the reader keeps in 

mind three basic "tricks" : we apply of the binomial theorem, 



(t-t') n = 5> fe (-t' 

k=0 



k( j.l\n—k 



(A7) 



we use Eq. (^) to simplify expressions, and we make changes of variables of the form 
x = $ (t'). These lead to expressions such as 



(d+1) 



'lnzo\ 



(n+l)/(d+l) 



/ 



t k (t') n ~ k dt' = [%{t)] k ^ d+ ^x {n - k ~ d)/{d+1) dx . 



(A8) 



Using Eq. fl3§ in Eq. (|A| 

d 



Vifa^) « 2tlv d (Jj {HI' 1 (t 2 - hf' 1 |m sti0 
Zo exp 



x 



t \ d+1 
K*b)(-j 



#1 1 |m st) i 

*o(*2) 



ljdt 

.d z o 

+ 



-Vb(ti,t 2 ) + ^^ I (lnz )- 1/(ci+1) (^ -hY-'ro 

fs>o(t 2 ) _ __, 

/ e dx 

2 (J) ItfrVst^-Wi,^) +^ rf ^ I (ln^ )- 1 /(' i + 1 )(t 2 -tO^Vo 
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Using Eq. in Eq. (gl 



$ (t 2 



7 



d + 1 



$o(*i) 



(A9) 



^o) 



TO 



x d 7 



d + 



d+1 



,*o(0 



2 |#rVst,o|f - $o(*) +^o(ln^o) d/( " +1) r - d 

' J 5TI.Vi)]-jf|( d i 1 )^-7 



x u 7 



d 



fc=0 
^0 



d + 1 
$o(*) 



*o(0 



dt' 



+ 



d + 



d-l 
- 1 fc=0 



d- 1 



fc/(d+i) 



.4 



d+1 

-(fc + 1) 1 
d + 1 'd+1 



,%(t) . (A10) 
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It is somewhat easier to evaluate $i(t). From Eq. fl42|) and Eq. 

= X(H)I Qvi f m (t')(t - t'Ydt' 
Jo 



wm (d+l) ln(zp) 
A (#) Id+i 



7" n 



+ / exp 
/o 



-ln(z )| 
A(F)Kt,o|<!-$o(t) + 



z (c?+ l)ln(2 ) 



X 



/;x;(f>*(-o-exp 



ln(zo) ( - 



dt' 



X(H)\m stt o\lz I 

\ U=o 

A(iy)|m sti0 | (z [E ( { ) [$o(t)] fe/(d+1) 7 (l - j^y, *o(f)) } 

- $o(*) J • (All) 



Combining Eq. Eq. ( [ATT] ), and Eq. (gSJ) gives mi(t) 
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FIG. 1. The relationship between the applied field H and system width L for a shorter (top 
solid curve) and a longer (bottom solid curve) fixed lifetime in a typical metastable magnetic 
system. Four regions are distinguished by differing decay processes: the Coexistence region (CE), 
the Single-Droplet region (SD), the Multi-Droplet region (MD), and the Strong-Field region (SF). 
The CE and SD regions are separated by the thermodynamic spinodal (dotted curve). The SD 
and MD regions are separated by the dynamic spinodal (dash-dotted curve). The SF region is 
separated from the other regions by the mean- field spinodal (dashed curve). [Reproduced from 



Fig. 1 of Ref. 33 
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FIG. 2. Configurations that may occur during the reversal process. As in the text, periodic 
boundary conditions are imposed, (a) A sketch of a "slab" configuration, (b) A typical realization 
of a single droplet in the process of overtaking the system. Grey squares are "up" spins and black 
squares are "down" spins. Here L = 60, # = -0.08 J, D = 0, T = 0.8T C , and £ = 410 MCSS. [Figure 
courtesy S. W. Sides.] (c) A typical realization showing the nucleation and growth of several 
droplets in the process of switching the magnetization. Here L = 120, H = —0.2J, D = 0, T = 0.8T C , 
and £ = 114 MCSS. [Reproduced from Fig. 5(b) of Ref. 133.] 
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FIG. 3. The approximate restricted free energy function F(m) as determined by Eq. (Q) with 
d = 2, T = 0.8T C , and L D = 500. (a) L = 5, H = 0.1J. (b) L = 10, H = 0.1J. (c) L = 500, H = 0. 
(d) L = 5000, H = 0. 
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FIG. 4. The relative lifetime vs. x = 2m sp D /\H\ in the SD region as given by Eq. (| 
|-ff| = 0.2J, T = 0.8T C , and L = 10. Each Monte Carlo point represents 47 500 decays. The quantity 
Hi was obtained in Ref. 33 from data at D = and is not the result of a fit to the -D-dependence. 
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FIG. 5. The metastable magnetization vs. H in the MD region. The data points were estimated 
by Ramos et a/.il by extrapolating back to t = assuming Avrami's Law. The smooth curve is a 
useful estimate, but the correct form of m ms vs. H is not known except near H = 0, where Eq. ( |1 lb| ) 
applies. The straight line indicates the approximation Eq. (|llb| ). 
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FIG. 6. The magnetization us. time in the MD region as given by Eq. ( |34| ) and Eq. (p8|). 
H = 0.3J, T = 0.8T C , and L = 50. The two values of D displayed are D = and D = 0.015. The solid 
curves are mo(t) + Dmi(t); the dotted curve (hardly distinguishable on the scale of this figure) is 
mr(t). Each Monte Carlo point represents 100 decays. 
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FIG. 7. The lifetime vs. D in the MD region for T = 0.8T C with H = 0.2J and # = 0.3 J. The 
solid curves represent the theoretical predictions given by combining Eq. (j5Q ) and Eq. (|5l| ) and 
are raoi the result of a fit to the D-dependence. Each Monte Carlo point represents at least 5 000 
decays in a system of size L = 100. [Reproduced from Fig. 4 of Ref. 34.] 
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